##############################################################
# REPLICATION CODE: "Are coups good for democracy?"          #      
# PURPOSE: Produce Figure 3 of the main text and Figure C-5  #
# of the Online Appendix, using repression_monthly.csv.      #      
# AUTHOR: George Derpanopoulos (gderpa@ucla.edu)             #
##############################################################

##############
## PREAMBLE ##
##############

# rm(list=ls())  # clear working space (recommended)
setwd('~/Desktop/R&P submission/replication files/repression_monthly')  # set working directory (adjust by user)

dat <- read.csv('repression_monthly.csv')  # load data

## Load packages (install any not installed already)
library(plyr)              # for manipulating data
library(ggplot2)           # for plotting
library(gridExtra)         # for multi-plot panels
library(stargazer)         # for LaTeX tables


#############
# FUNCTIONS #
#############

# For function below: breaks-down dataframe by cowcode and starts counting observations from 0, so if failed.coup occurs in first months of 1990, the subscript used in subsetting deaths.state becomes negative --> solution: drop 2 problematic cases
dat$failed.coup[dat$cowc=='NIG' & yearmonth=='1990-04'] <- 0
dat$failed.coup[dat$cowc=='ZAM' & yearmonth=='1990-06'] <- 0

## Function: for all coup events of chosen type: calculates within-country difference in total deaths by state-sanctioned perpetrators for 12 months before vs. 12 months after the coup; returns dataframe w. labels for plotting.
within.diff.fn <- function(coup.type,  # domain: {'failed.coup', 'dict2dict.coup', 'dict2dem.coup', 'reshuffling.coup'}
                           violence){  # domain: string of type of violence (e.g. 'deaths.state)
  within.diff <- ddply(dat, 
                       'cowc', 
                       function(x){
    wd <- numeric()
    coups <- which(x[, coup.type]==1)
    if(length(coups)>0){
      temp <- numeric(length(coups))
      if(coup.type=='failed.coup'){
        for(j in 1:length(coups)){
          if(x$democracy[coups[j]]==0){
            temp[j] <- sum(x[(coups[j]+1):(coups[j]+12),
                             violence], 
                           na.rm=T) -
              sum(x[(coups[j]-12):(coups[j]-1),
                    violence],
                  na.rm=T)
          }
          else{temp[j] <- NA}
        }
      }
      else{
        for(j in 1:length(coups)){
          temp[j] <- sum(x[(coups[j]+1):(coups[j]+12),
                           violence],
                         na.rm=T) -
            sum(x[(coups[j]-12):(coups[j]-1),
                  violence],
                na.rm=T)
        }
      }
      return(cbind(wd, temp))
      }
    else{}
  }
  )[,2]
  coup.type <- ifelse(coup.type=='failed.coup', 'F',
                      ifelse(coup.type=='dict2dict.coup', 'A',
                             ifelse(coup.type=='dict2dem.coup', 'D',
                                    ifelse(coup.type=='reshuffling.coup', 'N',))))
  return(data.frame(within.diff,
                    coup.type))
}



######################
# COUP TYPES/NUMBERS #
######################

# Dataframe w no. of coup events of each type; for plotting
coups <- data.frame(sum(dat$failed.coup==1 & dat$democracy==0, na.rm=T),
                    sum(dat$dict2dem.coup==1, na.rm=T),
                    sum(dat$dict2dict.coup==1, na.rm=T), 
                    sum(dat$reshuffling.coup==1, na.rm=T))
names(coups) <- c('Failed',
                  'Democratizing',
                  'Adverse regime change',
                  'No regime change')



######################
# FIG. 3 (MAIN TEXT) #
######################

# Dataframe of within-country changes in deaths by state for all coups/types; for plotting
within.diff.df <- data.frame(rbind(within.diff.fn(coup.type='reshuffling.coup', 
                                                  violence='deaths.state'),
                                   within.diff.fn(coup.type='dict2dict.coup', 
                                                  violence='deaths.state'),
                                   within.diff.fn(coup.type='dict2dem.coup', 
                                                  violence='deaths.state'),
                                   within.diff.fn(coup.type='failed.coup', 
                                                  violence='deaths.state')
))
drop <- which(within.diff.df$within.diff==4078)  # drop outlier to aid visualization
within.diff.df <- within.diff.df[-drop, ]  


# Figure 3 (main text)
pdf('within-deathsbystate.pdf')
ggplot(within.diff.df, 
       aes(x=coup.type,
           y=within.diff)) + 
  geom_boxplot(aes(fill=factor(coup.type)), alpha=.5) +
  geom_point(aes(color=factor(coup.type)), size=2.5) +
  labs(x="Coup type",
       y="Deaths(t) - deaths(t-1)") +
  theme(legend.position="none", text=element_text(size=20)) +
  coord_flip()
dev.off()



######################
# FIG C-5 (APPENDIX) #
######################

# Function for overlaying mean tseries in ggplot
stat_sum_single <- function(fun, color, geom='line') {
  stat_summary(fun.y=fun, color=color, geom=geom)
}

# Deaths by state for country-months around failed coups in dictatorships; for plotting
deaths.state.failed.dict <- subset(dat, 
                                   !is.na(months2.failed.coup)
                                   & democracy==0,  
                                   c(cowc, yearmonth, 
                                     months2.failed.coup, 
                                     deaths.state))

# ...jittered scatterplot w. mean line
deaths.state.failed.dict.plot <- ggplot(deaths.state.failed.dict, 
                                        aes(x=months2.failed.coup,
                                            y=deaths.state)) + 
  geom_jitter(color='red') +
  stat_sum_single(mean, geom='line', color='red') +
  geom_vline(xintercept=0) +
  labs(title=paste0('Failed', ' ', '(', as.numeric(coups['Failed']), ')'),
       x="Months since coup",
       y="Deaths")

# Deaths by state for country-months around dict2dict coup; for plotting
deaths.state.dict2dict <- subset(dat, 
                                 !is.na(months2.dict2dict.coup), 
                                 c(cowc, yearmonth, 
                                   months2.dict2dict.coup, 
                                   deaths.state))

# ...jittered scatterplot w. mean line
deaths.state.dict2dict.plot <- ggplot(deaths.state.dict2dict, 
                                      aes(x=months2.dict2dict.coup,
                                          y=deaths.state)) + 
  geom_jitter(color='purple') +
  stat_sum_single(mean, geom='line', color='purple') +
  geom_vline(xintercept=0) +
  labs(title=paste0('Adverse regime change', ' ', '(', as.numeric(coups['Adverse regime change']), ')'),
       x="Months since coup",
       y="Deaths")

# Deaths by state for country-months around dict2dem coup; for plotting
deaths.state.dict2dem <- subset(dat, 
                                !is.na(months2.dict2dem.coup), 
                                c(cowc, yearmonth, 
                                  months2.dict2dem.coup, 
                                  deaths.state))

# ...jittered scatterplot w. mean line
deaths.state.dict2dem.plot <- ggplot(deaths.state.dict2dem, 
                                     aes(x=months2.dict2dem.coup,
                                         y=deaths.state)) + 
  geom_point(color='green') +  # all 0, so jitter looks bad
  stat_sum_single(mean, geom='line', color='green') +
  geom_vline(xintercept=0) +
  labs(title=paste0('Democratizing', ' ', '(', as.numeric(coups['Democratizing']), ')'),
       x="Months since coup",
       y="Deaths")

# Deaths by state for country-months around reshuffling coup; for plotting
deaths.state.reshuffling <- subset(dat, 
                                   !is.na(months2.reshuffling.coup), 
                                   c(cowc, yearmonth, 
                                     months2.reshuffling.coup, 
                                     deaths.state))
drop <- which(deaths.state.reshuffling$deaths.state==4066)  # drop outlier to aid visualization
deaths.state.reshuffling <- deaths.state.reshuffling[-drop,]  

# ...jittered scatterplot w. mean line
deaths.state.reshuffling.plot <- ggplot(deaths.state.reshuffling, 
                                        aes(x=months2.reshuffling.coup,
                                            y=deaths.state)) + 
  geom_jitter(color='blue') +
  stat_sum_single(mean, geom='line', color='blue') +  # stat_sum doesn't work
  geom_vline(xintercept=0) +
  labs(title=paste0('No regime change', ' ', '(', as.numeric(coups['No regime change']), ')'),
       x="Months since coup",
       y="Deaths")

# Figure C-5
pdf('pooled-deathsbystate.pdf')
grid.arrange(deaths.state.failed.dict.plot, 
             deaths.state.dict2dict.plot, 
             deaths.state.dict2dem.plot, 
             deaths.state.reshuffling.plot,
             nrow=2, ncol=2, main='')
dev.off()



# ############################################
# # FIG 3 FOR WEIGHTED DEATHS (NOT REPORTED) #
# ############################################
# 
# # Dataframe of within-country changes in deaths by state weighed by no. incidents+campaigns from which they arise, for all coups/types; for plotting
# within.diff.w8.df <- data.frame(rbind(within.diff.fn(coup.type='reshuffling.coup',
#                                                      violence='w8.deaths.state'),
#                                       within.diff.fn(coup.type='dict2dict.coup',
#                                                      violence='w8.deaths.state'),
#                                       within.diff.fn(coup.type='dict2dem.coup',
#                                                      violence='w8.deaths.state'),
#                                       within.diff.fn(coup.type='failed.coup',
#                                                      violence='w8.deaths.state')
# ))
# drop <- which(within.diff.w8.df$within.diff==2045)  # drop outlier to aid visualization
# within.diff.w8.df <- within.diff.w8.df[-drop, ]  
# 
# # ...plot 
# ggplot(within.diff.w8.df, 
#        aes(x=coup.type,
#            y=within.diff)) + 
#   geom_boxplot(aes(fill=factor(coup.type)), alpha=.5) +
#   geom_point(aes(color=factor(coup.type)), size=2.5) +
#   labs(x="Coup type",
#        y="Deaths(t) - deaths(t-1)") +
#   theme(legend.position="none", text=element_text(size=20)) +
#   coord_flip()

##########################################################